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SUMMARY 


A procedure for the determination of airplane model structure from flight data 
is presented. 'The model is based on a polynomial spline representation of the aero- 
dynamic coefficients , and the procedure is implemented by use of a stepwise regres- 
sion. First, a form of the aerodynamic force and moment coefficients amenable to 
the utilization of splines is developed. Next, expressions for splines in one and 
two variables are introduced. Then the steps in the determination of an aerodynamic 
model structure and the estimation of parameters are discussed briefly. The focus 
of the paper is on the application to flight data of the techniques developed. Here, 
the parameters estimated from large -amplitude maneuvers are compared with a baseline 
set of parameter estimates from standard small-amplitude maneuvers and steady- 
state measurements. The model is further validated by comparing the predicted air- 
plane motion with actual measured time histories. It is thus shown that the proce- 
dure represents a further step toward the determination of a global model of an 
airplane from flight data. 


INTRODUCTION 

A procedure is outlined in reference 1 for the determination of airplane model 
structure from flight data, including nonlinear aerodynamic effects. This procedure 
focuses on finding the form of and parameters in aerodynamic model equations using 
a modified stepwise regression and several decision criteria. The aerodynamic func- 
tions are approximated by polynomials in airplane response and input variables. The 
procedure described was successfully applied to small-amplitude maneuvers. When 
applied to large-amplitude longitudinal maneuvers, the data were first partitioned 
into subsets as a function of angle of attack. Then, each subset was analyzed sepa- 
rately. This approach, however, has only limited application. First, each data 
subset must have a sufficient number of data points for successful model determina- 
tion. Second, for longitudinal maneuvers, the angle-of-attack intervals for individ- 
ual subsets must not be so small that the functions vary so little that accurate 
parameter estimation would not be possible. 

In large-amplitude and high-angle-of-attack maneuvers the behavior of aerody- 
namic functions in one region of angle of attack may be quite different from and 
totally unrelated to their behavior in another region. In these cases, the polyno- 
mial approximation for some aerodynamic nonlinearities would be inadequate. Polyno- 
mials are determined everywhere by their values in any interval, no matter how small. 
They can, therefore, follow a curve in one interval but depart from a curve or even 
oscillate widely elsewhere. Even if a higher-order polynomial approximates the aero- 
dynamic function sufficiently, the increase in the number of terms can lead to large 
covariances of their estimates. 

To avoid the disadvantages of the polynomial representation, spline functions 
can be used. Splines avoid some difficulties of polynomials because they are defined 
on preselected intervals and because the low-order terms may approximate various 
nonlinearities quite well. The application of splines to airplane model structure 
determination was first suggested in reference 2. Their use in real flight data 
analysis was then investigated, and the results are presented in references 3 to 6. 



The purpose of this investigation was to examine the use of polynomial splines 
in one and two variables for postulating the aerodynamic model equations and for 
determining a model structure by using a stepwise regression. This report is 
an extension of the research reported in references 1 , 5, and 6. The formulation of 
aerodynamic model equations and the definition of polynomial splines are discussed, 
followed by a discussion of model structure determination. The entire procedure is 
tested on several examples. 
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SYMBOLS AND ABBREVIATIONS 
coefficients of polynomial P(x 2 ) 
coefficient of ■ polynomial Q(x^) 
wing span, m 

general aerodynamic force and moment coefficient 

coefficient of x n in spline function of x 
h s 

coefficient of in s Pli- ne function of x.j,x 2 

rolling-moment coefficient, M^/qSb 
pitching-moment coefficient, M y /qSc 
yawing-moment coefficient, M^/qSb 

Ci 

longtiudinal-force coefficient, F x /qS 
lateral-force coefficient, F Y /qS 
vertical-force coefficient, F z /qS 
wing mean aerodynamic chord, m 

coefficient of (x - x„ . ) m in spline function of x. 

1 1i + « 

coefficient of (x - x - x . )** in spline function of x n ,x 9 

I 1 1 T Z ZT I + ' ^ 


ca qi 6e, l 


coefficients in spline function of C (a), C (a), and C (a), 


respectively 
F-statistic 


6e 


F , F ,F forces along longitudinal, lateral, and vertical body axes, respectively, N 
X Y Z 


g 

h 

k 

Z 


number of unknown parameters 
lag number 

number of spline knots over range of x 1 
number of spline knots over range of x 2 
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M maximum lag number 

M ,M ,M rolling, pitching, and yawing moments, respectively, N-m 
X Y Z 

m degree of polynomial spline in 

N number of data points 

n degree of polynomial spline in x 2 

P(x^) polynomial in x 2 of degree n 

p roll rate, rad/sec or deg/sec 

Q(x^) polynomial in x^ of degree m 

q pitch rate, rad/sec or deg/sec 

12 ,.. 

q = - P V , kinetic pressure, Pa 

2 

R squared multiple correlation coefficient 

r yaw rate, rad/sec or deg/sec 

2 

S wing area, m 

S (x) polynomial spline in x of degree m 

m 

S (x ,x ) polynomial spline in x 1 ,x 9 of degree m in x,. and degree n 
mn 12. 1 z 

m x 2 

2 

s estimated variance 

t time, sec 

V airspeed, m/sec 

W(h) autocorrelation function at lag h 

x general independent variable in one-variable spline approximation 

x(i) independent variable at time t. 


( Xj - x..) + - 


0 if x . < x. . 

1 

(x. - x. . ) m if x. > x. . 
„ D D 1 D li 


general dependent variable 
dependent variable at time t^ 
angle of attack, rad or deg 


midpoint of the a-interval of subset of partitioned data, rad or deg 
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0 

6 

a 

6 

e 

6 

r 

e(i ) 

0. 

J 

v ( i ) 

P 


angle of sideslip, rad or deg 
aileron deflection, rad or deg 
elevator deflection, rad or deg 
rudder deflection, rad or deg 
equation error at time t 

coefficient of general independent variable 
residual value at time t^ 

3 

air density, kg/m 


x. 


J 


a standard error of aerodynamic coefficient 

Superscripts: 

h,s degree of a polynomial 

• derivative with respect to time 

A estimate 

Subscripts: 
crit critical 

i or j value at which spline knots occur 
max maximum value 

p partial 

Abbreviations : 

ML maximum likelihood 

MSR modified stepwise regression 

Aerodynamic derivatives referenced to a system of body axes with the origin at 
the airplane center of gravity: 
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AERODYNAMIC MODEL EQUATIONS 


For the determination of model structure and estimation of aerodynamic param- 
eters, the analytical form of the aerodynamic model equations must be postulated. 

The expressions for coefficients of aerodynamic forces and moments used in this 
report are based on the following principal assumptions: 

1 . The instantaneous aerodynamic forces and moments depend only on the instanta- 
neous values of response and input variables. That is, no unsteady aerody- 
namic effects are considered. 


2. The dependence of longitudinal and lateral coefficients on response and input 
variables can be expressed as 


C a = 


(a = X, Z, or m) 


and 


c = C (p, 0 C,p,r ,6 ,6 ) 

a. G. Q L 


(a = Y, A, or n) 
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3, The resulting aerodynamic coefficients are obtained as sums of contributions 
due to (a, (3)/ and p, q, r, 6 a , 6 e , and S r . The second group of these 
contributions is, in general, a-dependent . 

Considering the preceding assumptions, the aerodynamic model equations can be 
written as 

c: = C (a, p)_ c _ n + C_ (a) qc/2V + C (a) 6 e (a = X, z, or m) (1) 

a a q-o e -u a q a §e e 


and 


C = C ( a, (3) n + C (a) pb/2V + C (a) rb/2V + C (a) 6 

a a K p=r=6 =6 =0 a a a c a 

a r p r oa 

+ C (a) 6 (a = Y, A, or n) (2) 

a c r 

6r 

The expressions for the aerodynamic coefficients are similar to those used in 
wind-tunnel testing practice. The first terms on the right-hand side of equa- 
tions (1) and (2) represent "static" parts with controls fixed at zero deflections. 
The remaining terms represent contributions of dynamic stability derivatives and 
control derivatives and their dependence on a. In equations (1) and (2), no 
a and (3 terms are explicitly introduced because of their near-linear dependence 
on the remaining variables. The effects of a and (3 are included primarily in 
contributions due to angular velocities. 

The form of equations (1) and (2) indicates that each term in these equations 
can be approximated by a spline either in (a,{3) variables or in the a variable. In 
longitudinal maneuvers with small lateral coupling, equation (1) can be further 
simplified by replacing the two-dimensional terms in (a,|3) by two terms in a. These 
equations then take the form 


C a = C a (a) P =q=6 e =0 + C a p 2 (a) P + C a q (a) ^/ 2V 


+ C *6e (a) 


(a = X, Z, or m) (3) 


Equations (2) and (3) represent fairly general formulation of aerodynamic coeffi- 
cients. In each particular case, however, the postulated aerodynamic model equations 
should reflect any available a priori knowledge, based on wind-tunnel and/or theoret- 
ical aerodynamic data. 


POLYNOMIAL SPLINES IN ONE AND TWO VARIABLES 

Spline functions are defined as piecewise polynomials of degree m. When con- 
tinuity restrictions are considered, the function values and derivatives agree at the 
points where the piecewise polynomials join. These points are called "knots" and are 
defined by the value of their projection onto the plane (or axis) of independent 
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variables. A polynomial spline of degree m with continuous derivatives up to 
degree m - 1 approximating a function f(x) for x G ^ X Q ,x max ^' ° an be ex P ressed 
as 


S (x) 
m 


m 


h=0 


C. x 
h 


k 

£ 

i=1 


D^x - x^) 


m 


where 


/ \ ^ 
(X - x.) + 



(x > x i ) 
(x < x ± ) 


(4) 


The values x 1 , . .., x fc are knots which obey the condition, 

x n < x, < x. < ... < x, < x , and CE and D, are constants. The spe- 
cial case of equation (4) for m = 0 (a spline of degree zero) represents 
an approximation by piecewise constants. 

The problem of multidimensional splines is addressed in reference 7. A space of 
these splines is constructed by taking the tensor product of one-dimensional spaces 
of polynomial splines. Because of the tensor nature of the resulting space, many of 
the simple algebraic properties of ordinary polynomial splines in one dimension are 
carried over. A spline in two variables x^ and x^ can be introduced for the 
approximation of a function f(x^,X 2 ) for x>| C t x io ,x 1max^ and x 2 ^ x 20 ' x 2max ^ • 
Then, as in the one-dimensional case, the two ranges t x -|o fX 1max^ and ^ X 20' X 2max^ are 
subdivided by sets of knots x^ ^ and x^^ where 


X 10 < x n 


< x „i < x „ 

Ik Imax 


20 


< x 


21 


< x n < x 

21 2max 


The points ( x 1 i^ x 2 i^ partition the above rectangle into rectangular panels. A poly- 
nomial spline of degree m for x^ and of degree n for X 2 with continuous 
partial derivatives up to (m - 1 ) + (n - 1 ) degree on the rectangle defined by the 

intervals [x 1Q ,x 1lnax ] and [ x 2 o' x 2max ] can be formulated as 


m n 


mn 


<*,'*2' - E. E c hs * £ VV <*, - *u>: + £ VV (x 2 

h=0 s=0 i=1 D=1 


.n 

- x„ . ) 
2j + 


e e v x < - «,!>; <* 2 -v+ 

1 = 1 J =1 


(5) 


where p ^( x 2^ and Qjtx^) are polynomials of degree n and m, respectively, and 


C hs and D ij are constants - 
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As examples of using splines in the approximation of aerodynamic functions, the 
vertical-force coefficient C z and yawing-moment coefficient C R are considered. 

In the first case, the form of C (a,q,6 ) based on splines can be, according to 
equation (3), written as 

C Z (a,q,6 e ) = C Z (a) q=S =0 + C z (a) qS/2V + C z (a) 6 e (6) 

e ‘ q °e q 6e e 


where 


Cl a) = c 0 + c a + Y D . a - ct. 
z z z ai i + 

a i=1 


C z - C z + Z D i<« - «i> + 

q q i=1 


(a) = C + Y D. . (a - a. )° 

Z c z _ A-* 6e,i x + 

6e oe x = 1 j 


(7) 


Equation (7) indicates that C z (a) is approximated by piecewise linear polynomials 
(the first-degree spline), whereas the remaining two functions are approximated by 
piecewise constants (the zero-degree splines). 

In the second case, 


c n («,p,p,r,6 a ,6 ir ) = c n (a,p) , 6 + C n (a) pb/2V 


a r 


+ C (a) rb/2V + C (a) 6 + C (a) 6 

n r n 6a a n 6r r 


( 8 ) 


Using equation (5) for x 1 = a and x 2 = f3, and selecting m = 0 and n = 1, the 
function C ( a, |3) can be approximated as 


k 

c n (a,p ) = C Q + + Y <A oi + A p)(a - a.) + 

i=1 


Z 


p. 


k l 


+ > B 0j T - mi + 

j=1 N / + 



D ij p 1 mi (a " V+ (9) 


i=i j=i 
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where 


<ipi < pj) 

(p > p.) 

(p < -Pj ) 

Pj . In this approximation of C n (a,f3), it is 
assumed that C n (p) is an odd function. The remaining functions in equation (8) are 
then represented by splines in a. 


MODEL STRUCTURE DETERMINATION 

The determination of an adequate model using the stepwise regression includes 
three steps: the postulation of terms which might enter the final model, the selec- 

tion of an adequate model, and the verification of the model selected. 

As shown in the previous section, the general form of aerodynamic model 
equations can be written as 


pi - 


ft* 


ei/ + 




There is always a positive value for 


y(t) = e Q + + 


,. + e x 

g-1 g-1 


( 10 ) 


In this equation, y(t) represents the resultant coefficient of aerodynamic force or 
moment (the dependent variable), 0^ to 0 are the constants in spline represen- 

tation of the aerodynamic functions, and x-j to Xg_i are the airplane response and 
input variables and their combinations (the independent variables). When the aero- 
dynamic model equations are postulated, the determination of significant terms among 
the candidate variables (determination of model structure) and estimation of corre- 
sponding parameters follow. 

Assuming that a sequence of N observations of y and of x at times 
t^, t^, *..., t N has been made, the measured values denoted by y(i) and x(i), 

where i = 1, 2, • • . , N, are related by the following set of N linear equations: 

y ( i) = 0^ + 0 x (i) + ... + 0 x (i) + e(i) (11) 

011 g-1 g-1 


where e(i) represents the equation error. An adequate model for the aerodynamic 
coefficients can be determined by applying the stepwise regression. The stepwise 
regression technique is described in references 1 and 8, and its main features are 
summarized in the appendix. 

When formulating spline functions in expressions for an aerodynamic coefficient, 
the degree of spline and the number and location of knots must be specified. A 
set of knots is fixed in advance. The number of candidate knots for each spline 


is limited only by the available computer memory. The stepwise regression proce- 
dure then selects only the knots which are associated with statistically signifi- 
cant parameters. In this way, a suboptimal number of knots and their locations 
are obtained. A more detailed discussion on the application of statistical variable 
selection technique to fit splines is contained in reference 9. The selection of 
degree of spline is influenced by the form of an aerodynamic function which is to be 
approximated. With little or no knowledge of that form, the model postulation can 
start with splines of a low degree, that is, the zero degree or first degree. After 
a tentative model structure is determined, the procedure can be repeated with higher- 
degree splines for better approximation of the data. 

Because of the combination of spline representation with the stepwise regression 
technique, the number of knots for each spline in the postulated model structure 
is limited only by practical considerations, that is, the available computer memory. 
The knots can be positioned arbitrarily within the range of independent variables. 

The estimation technique then selects only those knots associated with the important 
terms in the model equation. 

As explained in references 1 and 4, the regression analysis for model struc- 
ture determination and parameter estimation provides the opportunity to use subsets 
of measured data rather than the whole data set. This approach can, for example, 
simplify the analysis of maneuvers for obtaining the lateral parameters in equations 
for the coefficients C Y , C^, and C n . The measured data can be partitioned as a 
function of a, and it can be assumed that for each subset the coefficients are func- 
tions of a = , where is the midpoint of an a-interval of the ith subset. 

Then, for each subset, equation (2) can be simplified to 


C a (ct, |3,p,r,6 a ,6 r ) = C a (a,p) =6 =6 =0 + (a) pb/2V + (a) rb/2V + C (a) 6 a 

a r p r ua 

(a = Y, or n) (12) 


+ C (a) 6 
a 6r r 


Thus, the spline functions in two variables (a,p) are replaced by splines in p, and 
the remaining splines are replaced in a by constants. 

The last step in model structure determination and parameter estimation is model 
verification. The parameter estimates must have realistic values and should be com- 
pared with wind-tunnel results and theoretical predictions. Wherever possible, the 
least-squares estimates should be also compared with the estimates using different 
techniques, that is, the maximum-likelihood estimation method. Ultimately, the model 
should be a good predictor of airplane motion within the region of its assumed 
validity. 


EXAMPLES 

In the following examples the technique for model structure determination and 
parameter estimation was applied to measured data of a single-engine, low-wing 
research airplane. The basic characteristics of the airplane and instrumentation 
system are presented in reference 10. This airplane had undergone certain wing 
leading-edge modification which allowed the airplane to be trimmed at angles of 
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attack up to approximately 24°. The measured data were available in the form of 
input and output time histories sampled at 0.05-second intervals. The input vari- 
ables included the directly measured control-surf ace deflections. The variables a 
and p were measured by wind vanes, and the airspeed was measured by an anemometer. 
These measurements were corrected for the local airflow and offset with respect to 
the airplane center of gravity. The remaining output variables included angular 
rates and linear accelerations. These variables, together with the airspeed, were 
used to compute the aerodynamic coefficients (ref. 1), for which the angular acceler- 
ations were obtained by differentiation of spline fits of measured angular rates. 
These data included basically three different sets of maneuvers. For the first set, 
small-amplitude longitudinal and lateral maneuvers were excited by control-surface 
deflections at trimmed conditions within the range of a between 4° and 24°. From 
these maneuvers local models of aerodynamic coefficients were estimated. The second 
set included the data from longitudinal quasi-steady flights. These flights were 
represented by slow deceleration-acceleration maneuvers from which static parameters 
were determined. Finally, the third set of data consisted of large -amplitude longi- 
tudinal maneuvers and large -amplitude combined longitudinal and lateral maneuvers 
with the a variation between 0° and 30°. The large-amplitude maneuvers were ana- 
lyzed for determining the parameters of an extended model, that is, a model valid 
over an extended range of a. The combined maneuvers were intended for the determi- 
nation of the parameters of a model which would be valid within the flight envelope 
containing prestall and stall maneuvers. 

In figure 1, the three important longitudinal stability parameters estimated 
from 30 small-amplitude maneuvers are plotted against the angle of attack correspond- 
ing to the trimmed conditions. All these maneuvers were analyzed by the modified 
stepwise regression ( MSR) described in reference 1. Although the definition of "best 
model" is subjective and an exhaustive search of all candidate models is prohibitive 
in both cost and time, experience has shown that the MSR gives an adequate model. 

For the verification of results obtained, the maximum likelihood (ML) estimation 
technique of reference 11 was applied to nine test runs. In this analysis, the model 
structure was the same as determined by MSR. The ML estimates are represented in 
figure 1 by closed circles. Also plotted in figure 1 are the parameters obtained 
from quasi-steady flight using expressions from reference 10. 

The values of three lateral parameters obtained from 30 small-amplitude lateral 
maneuvers are given in figure 2. In this case the yawing-moment derivatives were 
selected. These parameters exhibit large scatter mainly in the stall and post-stall 
region. This scatter is caused by small excitation of yawing motion of the tested 
airplane. The estimates of the remaining lateral parameters are shown subsequently 
to be more consistent. As in the previous case, some runs were also analyzed by the 
ML method, and the results are presented in figure 2. The estimates of longitudinal 
and lateral parameters from small -amplitude maneuvers were assembled as a data base 
for the comparison with the results from large -amplitude maneuvers. 


Parameters From Large -Amplitude Longitudinal Maneuver 

A large-amplitude longitudinal maneuver was analyzed using the polynomial spline 

representation of aerodynamic force and moment coefficients. The time histories 

of the input and response variables in this maneuver are presented in figure 3. From 

these time histories the aerodynamic functions C - C^, and C were computed, and 

x ^ m 

they are plotted in figure 4 against a rather than time. In figure 5, the varia- 
tion of rate of pitch q and elevator deflection 6 0 with a is also shown. Both 
figures indicate the range of dependent and independent variables available for use 


in the regression equations and the distribution of measured points within the 
a-range. This information can be utilized for the selection of knots and degree of 
splines in the postulated model structure. 

For the spline approximation of all three aerodynamic coefficients, the form of 
equation (6) was used, first with the first-degree spline for C^fa), C z (a), and 
C m (a), and then with the zero-degree splines for the remaining terms. (See eqs. (7) 
for approximation of C z coefficient. ) Seventeen knots for each spline were pos- 
tulated as =6°, a 2 = 7 °' •••' a i 7 = 22 ° # The resulting estimates showed that 

the zero-degree spline approximation of q- terms was rather coarse. Therefore, the 
first-degree and second-degree splines were introduced instead. The second -degree 
spline was then considered as the final approximation. As an example of options used 
in the q-terms representation, the C (a) estimates are presented in figure 6. The 

m q 

different spline approximations of all q-terms had only a small effect on the fit to 
the data and on the estimated parameters in the remaining terms. The final estimates 
of polynomial spline terms representing the three aerodynamic coefficients are plot- 
ted against a in figure 7. The aerodynamic function estimates are compared with 
the estimates from small-amplitude maneuvers and show very good agreement with those 

results. The fit of the polynomial and spline models to the data and the usefulness 

of the terms in those models can be revealed by the comparison of standard error 
of aerodynamic coefficient a and the squared multiple correlation coefficient R 2 . 
(See appendix for more detailed explanation.) The values of both coefficients are 
summarized in table I. For small-amplitude maneuvers, the low values of a and high 

values of R 2 correspond to low-angle-of -attack regimes? conversely, a is larger 

and R 2 smaller for values of a between 20° and 24° • 

In figure 8, the spline terms of C z (a) and functions are compared with 

the measurement of these relationships in the quasi-steady flight. This measurement 

resulted in a double-value function C„ ( a) and C (a) for values of a between 10° 

Zj 4 in 

and 22°, depending on increasing or decreasing values of a. This phenomenon can 
be caused by the aerodynamic hysteresis and by the hysteresis in the control system. 
Because of the relatively small differences in both branches of C z (a) and C m (a) 
curves with respect to the accuracy of these estimates, the hysteresis was not 
modeled in equation (6). 

Even if the agreement between the results from small- and large-amplitude maneu- 
vers is very good, the resulting model is further verified by simulating the airplane 
longitudinal response. This is done by using the extended model approximated by 
splines and by using the elevator deflection time history from a selected independent 
maneuver. The time histories of input and response variables are presented in fig- 
ure 9, and the response variables V, a, and q are compared with those predicted 
by the model. The comparison reveals good prediction capabilities for the model. 


Parameters From Large -Amplitude Combined Maneuvers 

The time histories of one of the combined maneuvers are given in figure 10. The 
response variables exhibit a persistently excited lateral motion due to the rudder 
and aileron and longitudinal motion due to elevator deflections and coupling between 
lateral and longitudinal modes. From this particular maneuver, the variation of 
aerodynamic coefficients with a is shown in figure 1 1 , and the variation of lateral 
variables is shown in figure 12. Both figures indicate the amount of excitation of 
dependent and independent variables in regression equations. 
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In the first step of data analysis, simplified models for the lateral coeffi- 
cients were formulated. For the coefficient C , the aerodynamic model equation had 
the form 


C 

n 


C ( a) 


p + C (a) pb/2V 

n p 


+ C (a) rb/2V + C 


( a) 


6a 


6 + C (a) 6 

a n P r 

6r 


(13) 


with zero-degree splines for all terms included. The equations for C Y and had 

a similar form. In the next step of the analysis, more complicated models with two- 
variable splines in a and p were used. This formulation is expressed by equa- 
tions (B) and (9) for the coefficient C • The main differences in the results of 

n 

the two approaches were found in parameters of the yawing-moment equation. This 
might indicate that the most pronounced nonlinearities in p are in the coefficient 
. Figure 13 presents the estimates of 


n 


p 


within the range of p from -4° 


to 4° from the three combined maneuvers. Each maneuver covered a different range of 
a and used a C_ (a) and C (a,p) spline in the model. Both sets of estimates are 


n p 


n 


compared with those obtained from small-amplitude maneuvers. As shown in figure 13, 
the second approach gives more consistent estimates which also agree better with the 
data base than the results of the first approach. The only difference remains in the 
region of a from 16° to 24°. In this area, the spline approximation of c n ( a 'P) 
resulted in directional instability, whereas the small-amplitude data have essen- 
tially positive values. Exclusion of the two-dimensional spline c n (a,p) from the 
model also had an effect on the remaining parameters, mainly the damping in yaw C n • 


The estimates of this parameter are given in figure 14. The second model improved 
the consistency of the estimates and the agreement with the results of small- 
amplitude manuevers . 


The estimates of three parameters in the rolling-moment equation with the spline 
C^(a, p) in the model are shown in figure 15. Consistent results were obtained 
for parameters C and C , but some scatter is apparent in the estimates of C 

l r 


From the remaining lateral parameters 

were estimated with good confidence, 
estimates of the other parameters. 


, the values of , 


C , and C 
h-' n 6r 


P or 6a 

but rather larger scatter was observed in the 


For the verification of previous results and for obtaining more accurate esti- 
mates of lateral parameters, the measurements from 12 combined maneuvers were joined 
together into one set of data. The resulting ensemble of about 13 000 data points 
was then partitioned into 22 subsets according to the values of a. The modeling 
of the lateral parameters was conducted mostly on 1° subspaces of the 0° to 30° 
a-space. As an example, the model for was postulated as 


C n (-a,p, P ,r,6 a ,6 r ) = C n P + Y] V. 

P i^1 P 1 


P 1 " 


Pi 


PI A n 


+ C pb/2V 


+ C n rb/2V + C n 6, + C n 6_ 
n r n 6a a n 6r r 


(14) 
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where a denotes the midpoint of an a-interval. The knots of the spline in |3 were 
selected at 4°, 8°, 12°, 1 6° , and 20°. The estimates of parameters from partitioned 
data are presented in figure 16. In general, these estimates are more consis- 
tent than those from individual maneuvers and are closer to the results from small- 
amplitude maneuvers. The nonlinear variation of the lateral coefficients with |3 
for a = a, and with the remaining lateral variables equal to zero, is demonstrated 
in figure 17 for three different values of a. Finally, the values of C n (a,f3) for 
P > 0 and all remaining independent variables equal to zero are plotted in fig- 
ure 18. The results from the partitioned data have confirmed the previous conclu- 
sion about the most pronounced nonlinearities in with (3 and about the small 

nonlinear variation of the remaining two coefficients with this variable. 

For the additional comparison of results from small- and large -amplitude maneu- 

A 

vers and from individual runs and partitioned data, the values of standard error a 
and squared multiple correlation coefficient R 2 for all cases are summarized in 
table II. For small-amplitude maneuvers and partitioned data the minimum values of 

A o 

a and maximum value of R correspond to low values of a, and the other values 
correspond to high values of a. From the values of a and R 2 , it could be con- 
cluded that at high angles of attack only the model for was not able to fully 

explain all the variation in the measured data. This conclusion indicates that some 

additional terms in the model should be considered. Models for C used in large- 

n o .a 

amplitude maneuvers without partitioning brought some improvements in R and a 
over the same values from previous maneuvers. 

The combined maneuvers were primarily intended for the estimation of lateral 
parameters. However, an attempt was also made to estimate the parameters in the 
pitching-moment equation from these data. A possibility of aerodynamic coupling was 
considered and the model was postulated as 


C m (a,p,q,r,6 a ,6 e ,6 r ) 


= C + C qc/2V + C 6 

m p=q=r=6 a =6 e =6 r =0 m q m 6e e 

3 4 

+ c 6 + c |p| + £ c ip - p. | 

m c 3 e m „ 1 v ' , m . I r "l I + 

6e P i=1 pi 


+ C | pb/2V[ + C | rb/2V| + C |6 | + C 


m 


m c i a I 
6a 


m c I r I 
6r 


(15) 


where 


r 

o 


- PiU-o- 


p + p ± 


(I PI < p ± > 

(p > p.) 

(p < -pu 


The estimates of the most significant parameters in the model are presented in fig- 
ure 19. Their comparisons with the estimates from small- and large -amplitude longi- 
tudinal maneuvers indicate good estimates of C and C terms. Some discrepan- 

m m 6e 


14 


For the values of 


cies with previous results are visible in the estimates of C 

m q 

2 

C m from partitioned data, R varied between 63 percent and 37 percent. This indi- 
cates a strong possibility of modeling error in equation (15). The low values of R 2 
are in contrast with values of R 2 in table I for both small- and large-amplitude 
longitudinal maneuvers. 

For the verification of the estimated lateral parameters, the aerodynamic model 
from partitioned data was used. The equations of motion for a = 22.5° were inte- 
grated with the initial conditions and control time histories from a flight trimmed 
at a « 20° • The control input for this independent maneuver consisted of a doublet 
in aileron followed by a doublet in rudder. The predicted time histories are plotted 
against the actual flight data in figure 20. The comparison of these time histories 
shows an acceptable agreement between them. 


CONCLUDING REMARKS 

The application of polynomial splines and stepwise regression for the determina- 
tion of airplane model structure from flight data has been demonstrated. First, a 
form of the aerodynamic force and moment coefficients compatible with the utilization 
of splines in one (angle of attack) and two (angles of attack and sideslip) vari- 
ables was developed. Then, the model postulated for the analysis of flight data from 
large-amplitude maneuvers was discussed. The procedure of model determination was 
used in several examples. Whenever possible, the resulting parameters were verified 
by comparison with a baseline formed by the estimates from small-amplitude maneuvers 
and steady-state measurements. The airplane equations of motion, employing the model 
obtained, were numerically integrated to simulate the airplane responses to actual 
flight inputs. These responses were then compared with the measured flight time 
histories of response to that input. 

The main conclusions from the research reported can be summarized as follows: 

1 . In order to obtain a global aerodynamic model from flight data, it was advan- 
tageous to use large-amplitude longitudinal and combined maneuvers in which the lon- 
gitudinal and lateral motion were excited over an extended range of angle of attack. 

2. When postulating the aerodynamic model with little or no a priori knowledge 
of its form, the approximating splines can be of low degree (zero or first) . After 

a tentative model structure is determined, the procedure can be repeated with higher- 
degree splines for better approximation of the data. 

3. Because of the combination of spline representation with the stepwise regres- 
sion, the knots can be positioned arbitrarily within the range of independent vari- 
ables. The number of candidate knots for each spline is limited only by the avail- 
able computer memory. 

4. A set of one-dimensional splines in the angle of attack can closely approxi- 
mate uncoupled longitudinal maneuvers. For the combined maneuvers, however, two- 
dimensional splines in the angle of attack and sideslip must, in general, be con- 
sidered. To avoid the complexity of two-dimensional splines, the data from combined 
maneuvers can be partitioned into subsets according to the values of angle of attack. 
The lateral parameters from partitioned data were more consistent and closer to the 
baseline than results from individual large maneuvers. 


15 



The procedure presented represents another step towards the determination of 
a global model of an airplane from flight data. It can provide valuable informa- 
tion about aerodynamic forces and moments acting upon the airplane in large -amplitude 
maneuvers and/or maneuvers in high-angle -of -attack regions. 


Langley Research Center 

National Aeronautics and Space Administration 
Hampton, VA 23665 
January 17, 1983 



APPENDIX 


STEPWISE REGRESSION 

This appendix describes the basic principles and features of the stepwise 
regression which is used to determine aerodynamic model structure from flight data. 
This procedure begins with the assumption that there are no variables in the postu- 
lated regression equation other than the bias term 0 Q . An effort is then made to 
find an optimal subset of variables by inserting independent variables into the model 
one at a time. The first independent variable selected for entry into the equation 
is the one that has the largest correlation with the dependent variable y. Suppose 
that this variable is x^ . This is also the variable that produces the largest 
value of the F-statistic for testing the significance of regression. The variable is 
entered if the partial F-statistic exceeds a preselected critical F-value 


F 

P 


e 


i 


2 A 

8 ^( 0 ^ 


> F 


crit 


/V ^ A 

where 0-] is the estimated parameter associated with x^ and s ( 0 ^ ) is the vari- 
ance estimate of 0^ . 

The second variable chosen for entry is the one that now has the largest corre- 
lation with y after adjusting for the effect on y of the first variable entered, 
x^ in this case. These correlations are referred to as partial correlations. In 
general, at each step, the independent variable having the highest partial correla- 
tion with y is added to the model if its partial F-statistic exceeds the prese- 
lected F . , . At each step of the procedure, all variables entered into the model 
crit 

previously are also reassessed by examining their partial F-s tatis tics . A variable 
added at an earlier step may be redundant because the relationship between it and the 
remaining variables now in the equation has reduced its value of Fp to less than 
F crit* If th is happens, the significant variable is deleted from the regression 
model. The procedure terminates when all significant terms have been included in the 
model. 


As a new variable enters the model, several useful quantities are calculated at 
each stage of the stepwise regression. All these quantities should be examined for 
the final model selection. First, the user can consider the total F-value for a 
given model of Q variables calculated as the ratio of the mean square due to the 
regression to the mean square of the residual. This ratio is given as 

N 

23 iy(i) - yi 2 
p _ i= i N - g 

F N g - 1 

23 ty(i) - y(i) ] 2 

i=1 

where 

1 N 

y=i I y (i) 

i=1 
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This number usually increases to some maximum value as new variables enter the 
regression, but then decreases slightly as the new terms are less effective in 
reducing the residuals, Heuristically , the maximum F-value represents a model which 
best fits the data with a minimum number of parameters. Second, the squared multiple 
correlation coefficient R 2 is calculated. This number, expressed as a percentage, 
is a measure of the usefulness of the terms, other than 0 Q , in the model. The value 
of R 2 would be 100 percent for a model that perfectly fit the data. Third, at each 
stage, the partial F-values Fp for each parameter are printed. The user should 
look for consistency in the values of Fp. For example, if one value of is only 

slightly greater than F cr it and all other values of F^ are much greater, the user 
may not want to include the variable with the small value of F p in the model. The 
fourth aid in model selection is the estimated normalized autocorrelation function 
for the residuals. The estimate of the autocorrelation function at lag h is given 
by 


N-h 

W(h) ~ x. 1 u H v(i ) v(i + h) 
N-h ' 
i = 1 


(h = 0, 1 , • • • , M) 


where h is the lag number and M is the maximum lag number, which is usually 
10 percent of N. The normalized autocorrelation function is calculated as 
W(h)/W(0). This function should approach that for white noise with a value of 1 at 
zero lag and values of 0 at lag for 1 to M. In applications, when the value of F 
for a parameter makes the utility of an independent variable questionable, the con- 
tribution of that variable to the actual model structure can be assessed by observing 
the effect of the variable on the autocorrelation function of residuals. The fifth 
number that aids the user is the standard error in the residuals a, which is printed 
at each stage of the regression. 


One learns from experience that not all of the five criteria listed above are 
"optimally" satisfied for any single model. However, the stepwise regression and its 
associated information criteria do significantly reduce the number of possible models 
from which the user must choose. Moreover, as the model structure is determined, 
so are the parameter estimates. Finally, ambiguity in the model selection can also 
be resolved by requiring that the estimated parameters make sense physically and that 
the selected model have good prediction capability. 

The selection of a set of candidate model variables from which the stepwise 
regression can build a model should rely on the user’s a priori knowledge of the 
physical system that is to be modeled. For the airplane, such assumptions as the 
most influential variables and symmetry considerations have led to the following 
logic for selection of candidate model variables for a spline analysis of the longi- 
tudinal maneuver. (See eqs. (6) and (7) for C . ) Though it appears lengthy and 

z 

awkward, this formulation of the FORTRAN code allows for simple deletion, addition, 
and/or change in candidate model variables. 


DO 910 1=1 ,NPTS 
X(1 , I ) =ALPH( I ) 

X(2,I)=C/(2*VEL(I) )*Q(I) 

X( 3 , I=DELE( I) 

DO 91 1 111=4,39 
911 X(III,I)=0. 

IF( ALPH( I) .GE.XKN0T( 1 ) ) X( 4 , I ) =ALPH( I ) -XKN0T( 1 ) 
IF(ALPH(I) .GE.XKNOT( 1 ) ) X( 5, I) =X( 2,1) 
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IF ( ALPH( I ) ,GE.XKN0T(2) ) 
IF( ALPH( I) .GE.XKNOT( 2) ) 
IF(ALPH( I) .GE.XKNOT(3) ) 
IF(ALPH(I) ,GE.XKNOT(3) ) 
IF( ALPH( I) ,GE.XKNOT(4) ) 
IF(ALPH(I) ,GE.XKNOT( 4) ) 
IF(ALPH( I) ,GE.XKNOT( 5) ) 
IF(ALPH(I) .GE.XKN0T(5) ) 
IF(ALPH(I) .GE.XKN0T(6) ) 
IF(ALPH(I) .GE.XKN0T(6) ) 
IF(ALPH(I) ,GE.XKN0T(7) ) 
IF(ALPH(I) . GE . XKNOT ( 7 ) ) 
IF(ALPH(I) .GE.XKNOT(8) ) 
IF(ALPH(I) .GE.XKNOT(8) ) 
IF(ALPH(I) .GE.XKNOT(9) ) 
IF(ALPH( I) .GE.XKNOT(9) ) 
IF(ALPH(I) .GE.XKNOT(IO) 
IF ( ALPH( I) .GE . XKNOT ( 1 0) 
IF(ALPH(I) .GE .XKNOT ( 1 1 ) 
IF ( ALPH( I) .GE . XKNOT ( 1 1 ) 
IF(ALPH(I) .GE.XKNOT( 1 2) 
IF(ALPH(I) .GE .XKNOT ( 1 2) 
IF( ALPH( I) ,GE.XKNOT( 1 3) 
IF(ALPH(I) .GE .XKNOT ( 1 3) 
IF(ALPH(I) ,GE.XKNOT( 14) 
IF(ALPH(I) .GE .XKNOT ( 1 4) 
IF(ALPH(I) .GE .XKNOT ( 15) 
IF ( ALPH( I) .GE . XKNOT ( 1 5) 
IF(ALPH(I) .GE ,XKNOT( 16) 
IF(ALPH(I) .GE .XKNOT ( 1 6) 
IF(ALPH(I) .GE.XKNOT(1 7) 
IF(ALPH(I) .GE . XKNOT ( 17) 
IF(ALPH(I) ,GE.XKNOT(7) ) 
IF(ALPH(I) ,GE .XKNOT ( 1 3) 
910 CONTINUE 


X( 6 , 1 ) =.ALPH d ) -XKNOT ( 2 ) 

X( 7 , 1) =X( 2,1) 

X( 8 , 1 ) =ALPH( I) — XKNOT( 3 ) 

X(9 , I)=X( 2,1) 

X( 1 0 , 1 ) =ALPH( I ) — XKNOT( 4) 

X(1 1 ,I)=X(2,I) 

X( 1 2 , 1 ) =ALPH ( I ) — XKNOT( 5 ) 

X( 1 3 , 1 ) =X( 2,1) 

X( 14,l)=ALPH(I)-XKNOT(6) 

X( 1 5, 1) =X( 2,1) 

X( 1 6, 1 ) =ALPH( I ) -XKNOT ( 7 ) 

X( 1 7, 1) =X( 2,1) 
X(18,I)=ALPH(I)-XKNOT(8) 

X( 1 9 , 1) =X( 2,1) 

X( 20, 1 ) =ALPH( I) — XKNOT( 9) 

X( 21 , 1 ) =X (2,1) 

) X( 22, 1 ) =ALPH( I) — XKNOT( 10) 
) X( 23 , I) =X( 2,1) 

) X( 24, 1 ) =ALPH( I) -XKNOT( 11) 
) X( 25 , I) =X( 2,1) 

) X( 26 , 1 ) =ALPH( I ) — XKNOT ( 12) 
) X( 27 , I ) =X( 2,1) 

) X( 28, 1 ) =ALPH( I ) -XKNOT ( 13) 
) X( 29 , I) =X( 2,1) 

) X( 30,I)=ALPH(I)-XKNOT(1 4) 
) X( 31 , 1 ) =X( 2,1) 

) X( 32,I)=ALPH(I)-XKN0T( 1 5) 
) X( 33 , I) =X( 2,1) 

) X( 34,1) =ALPH( I) -XKNOT ( 1 6) 
) X( 35 , I) =X( 2,1) 

) X( 36,I)=ALPH(I)-XKNOT(1 7) 
) X( 37 , I ) =X( 2,1) 

X( 38 , I ) =X( 3,1) 

) X( 39, I)=X( 3,1) 


In the preceding printout, VEL(I) = Airspeed V at t ± , Q = Pitch rate q, 

NPTS = Number of data points N, C = Wing mean aerodynamic chord c, and 
X(J,I) = Value of jth model variable at d . The symbols XKN0T( ) indicate knots 
for specific values of a. The table actually gives the logic for creating the 
(39 x N) matrix containing the time histories of each of the 39 candidate independent 
variables. The 17 knots in angle of attack can be set at any value the user deems 
adequate for the data by setting XKNOT(I) in the program with I = 1,17. Changing 
the candidate model variables can easily be accomplished by substituting the new 
variable for any of the 39 candidates listed. The number of candidate variables is 
limited only by the size of the computer memory. 
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TABLE I.- COMPARISON OF STANDARD ERRORS AND SQUARED MULTIPLE CORRELATION 
COEFFICIENTS FROM DIFFERENT MANEUVERS FOR LONGITUDINAL DATA 




Small -amplitude 


Large -amp 1 i tude 

Aerodynamic 


maneuvers 


maneuvers 







coefficient 

a 

R 

2 

A 

o 






a 

R 2 


Min. 

Max. 

Min. 

Max. 



C X 

0.0011 

0.020 

71 .5 

99.9 

0.014 

96.5 

c z 

.0064 

.042 

88.3 

99.9 

.039 

98.6 

C m 

.0052 

.028 

89.1 

99.5 

.033 

93.5 


TABLE II,- COMPARISON OF STANDARD ERRORS AND SQUARED MULTIPLE CORRELATION 
COEFFICIENTS FROM DIFFERENT MANEUVERS FOR LATERAL DATA 



Small-amplitude 



Large-amplitude 





maneuvers 





maneuvers 




Aerodynamic 





Separate 

runs 


Partitioned data 

coefficient 

a 






• 


■ 








_ 

2 

- 


2 







y 

R 


y 

R 


Min. 

Max. 

Min. 

Max. 

Min. 

Max. 

Min. 

Max. 

Min. 

Max. 

Min. 

Max . 


0.0026 

0.0120 

94.3 

99.9 

0.0072 

0.0390 

96.4 

99.1 

0.0066 

0.0360 

81 .6 

99.3 

C, 

.0014 

.0062 

88.9 

96.5 

.0049 

.0084 

93.7 

94.5 

.0026 

.0095 

86.1 

97.4 

C n 

.0008 

.0040 

60.6 

99.3 

.0032 

.0057 

88.9 

96.4 

.0027 

.0064 

78.5 

96.6 
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O MSR estimates 


• ML estimates 



a, deg 


Figure 2.- Estimated lateral parameters from small-amplitude maneuvers. 

Transient data. 
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0 5 10 15 20 25 30 


a , deg 

Figure 4.- Time histories of longitudinal aerodynamic 
coefficients in large-amplitude maneuver plotted 
against angle of attack. 
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a , deg 


Figure 5.- Variation of two longitudinal variables with angle of attack 

in large-amplitude maneuver. 













Figure 8.- Comparison of vertical-force and pi tching-moment coefficient in steady 
state conditions estimated from different maneuvers. 



50 


MEASURED 
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Time, sec 

Figure 9.- Time histories of measured longitudinal flight data and 

those computed by using parameters obtained by stepwise regression 
with spline approximation. 
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Figure 10.- Time histories of measured input and response variables in large- 

amp litude maneuver. 
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Figure 12.- Variation of lateral variables with angle 
of attack in large -amplitude maneuver. 
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Figure 12.- Concluded. 
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Figure 14#- Comparison of damping-in-yaw derivative estimated from different 
maneuvers and using different spline approximations. 
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Figure 16.- Concluded. 


II III II ■Mill llll 







Ill■ll■llll■■l■l III II I 


ill i 111 hum ii iiHii ii iiiiii mi mu hi mi i mi mi mu n 


Will II III 1 1 ii mi I II I 


I I I 


nmiiHiiiiiii i n m ii iiiiii 



Figure 


18 .- 


Steady-state values 
amplitude 


of yawing-moment coefficients estimated from large- 
maneuvers. Partitioned data. 
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Figure 20,- Time histories of measured lateral flight data and those computed by 
using parameters obtained by stepwise regression from partitioned data. 
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